library(tidyverse, verbose = FALSE)
Registered S3 methods overwritten by 'dbplyr':
  method         from
  print.tbl_lazy     
  print.tbl_sql      
-- Attaching packages --------------------------------------------------------------------------------------------------------------- tidyverse 1.3.1 --
v ggplot2 3.3.5     v purrr   0.3.4
v tibble  3.1.3     v dplyr   1.0.7
v tidyr   1.1.3     v stringr 1.4.0
v readr   2.0.1     v forcats 0.5.1
-- Conflicts ------------------------------------------------------------------------------------------------------------------ tidyverse_conflicts() --
x dplyr::filter() masks stats::filter()
x dplyr::lag()    masks stats::lag()
FUN.read <- function(path, n){
        M <- numeric()
        for(i in 1:n){
                spec <- read_tsv(file = paste(path, "a",i,".ols",  sep = ''), 
                                 skip = 6, show_col_types = FALSE) %>% 
                        select(Counts)
                M <- rbind(M, spec$Counts) 
        }
        M
}

M_espectros <- vector(mode = 'list', length = 3)

dir <- "./espectros/cupla Zr-Er-ZrNb paralela/ExpI/"
n_spec <- length(list.files(dir))
M_espectros[[1]] <- FUN.read(dir, n_spec)

dir <- "./espectros/cupla Zr-Er-ZrNb paralela/ExpII/"
n_spec <- length(list.files(dir))
M_espectros[[2]] <- FUN.read(dir, n_spec)

dir <- "./espectros/cupla Zr-Er-ZrNb paralela/ExpIII/"
n_spec <- length(list.files(dir))
M_espectros[[3]] <- FUN.read(dir, n_spec)

# cantidad de archivos por barrido
lapply(M_espectros, dim)
[[1]]
[1]  235 7865

[[2]]
[1]  234 7865

[[3]]
[1]  238 7865
# wavelength values
wavelen <- read_tsv(file = paste(dir, '/a1.ols', sep = ''), skip = 6, show_col_types = FALSE) %>% 
                select(Wavelength) %>% 
                rowid_to_column()

wavelen$Wavelength <- round(wavelen$Wavelength, 4)
# La normalizacion es por detector y se realiza por suma total. Los indices 
# correspondientes a cada detector son:
#         
#       - 1:2048
#       - 2049:3983
#       - 3984:5924
#       - 5925:7865

# Funcion de normalizacion
FUN.norm <- function(M, n = 4, index){
        lista <- vector(mode = 'list', length = n)
        for(i in 1:n){
                # separar detector
                lista[[i]] <- M[,index[[i]]]
                # hacer minimo = 0
                minimo <- apply(lista[[i]], 1, min)
                lista[[i]] <- lista[[i]] + abs(minimo)
                # Dividir por suma total
                total <- apply(lista[[i]], 1, sum)
                lista[[i]] <- lista[[i]] / total
        }
        lista       
}

num_detec <- 4    # Numero de detectores
index_detec <- list(c(1:2048), c(2049:3983), c(3984:5924), c(5925:7865)) 

M_norm <- M_espectros %>% map(FUN.norm, n = num_detec, index = index_detec)

# suma el el espectro que fue separado para normalizar por detector
M_norm <- M_norm %>% map(~ cbind(.x[[1]], .x[[2]], .x[[3]], .x[[4]]))

lapply(M_norm, dim) # para inspeccion
[[1]]
[1]  235 7865

[[2]]
[1]  234 7865

[[3]]
[1]  238 7865
# Funcion que integra el pico
FUN.cuentas <- function(p, barrido, spec_ini = 1, spec_fin = 200){
        ind <- wavelen$rowid[which(wavelen$Wavelength == p)]
        m <- M_norm[[barrido]][spec_ini:spec_fin,(ind-2):(ind+2)]
        cuentas <- apply(m, 1, sum)
        cuentas
}

# Funcion para graficar perfil
# 
FUN.perfil <- function(df){
        g <- df %>% 
                pivot_longer(Er:Nb,
                             names_to = "Elementos",
                             values_to = "Intensidad") %>% 
                ggplot(aes(x = rowid*40, y = Intensidad, colour = Elementos)) +
                geom_point() #+ geom_line()
        
        caption <- paste("Lineas elegidas:","\n",
                         "Er: ",pico_Er,"\n", 
                         "Zr: ",pico_Zr,"\n", 
                         "Nb: ",pico_Nb,"\n", sep = '')
        
        g + annotate(geom = 'text', x = 20, y = 0.009, label = caption) +
                xlab('Distancia (micrones)') +
                ylab('Cantidad de cuentas (UA)')
}

pico_Er <- 369.0904
pico_Nb <- 322.3791
pico_Zr <- 357.0441
###### ## ##  PERFIL BARRIDO 1   ## ## ######
b <- 1
ini <- 19   #descarta los primeros 18 perfiles
fin <- 234  #toma hasta perfil 234
Er <- FUN.cuentas(pico_Er, barrido = b, ini , fin )
Zr <- FUN.cuentas(pico_Zr, barrido = b, ini , fin )
Nb <- FUN.cuentas(pico_Nb, barrido = b, ini , fin )

df_1 <- tibble(Er = Er, Zr = Zr, Nb =Nb) %>% rowid_to_column()
g_barrido_1 <- FUN.perfil(df_1)
g_barrido_1

b <- 2
ini <- 21
fin <- nrow(M_norm[[b]])
Er <- FUN.cuentas(pico_Er, barrido = b, ini , fin )
Zr <- FUN.cuentas(pico_Zr, barrido = b, ini , fin )
Nb <- FUN.cuentas(pico_Nb, barrido = b, ini , fin )

df_2 <- tibble(Er = Er, Zr = Zr, Nb =Nb) %>% rowid_to_column()
g_barrido_2 <- FUN.perfil(df_2)
g_barrido_2

b <- 3
ini <- 22
fin <- nrow(M_norm[[b]])
Er <- FUN.cuentas(pico_Er, barrido = b, ini , fin )
Zr <- FUN.cuentas(pico_Zr, barrido = b, ini , fin )
Nb <- FUN.cuentas(pico_Nb, barrido = b, ini , fin )

df_3 <- tibble(Er = Er, Zr = Zr, Nb =Nb) %>% rowid_to_column()
g_barrido_3 <- FUN.perfil(df_3)
g_barrido_3

¿Cual será la referencia para unir los datos?

# desplazando indices
# rowid viene a ser el giro del micrometro. 40 micrones
df_2$rowid <- df_2$rowid + 3
df_3$rowid <- df_3$rowid + 6

# esto hay que cambiarlo a distancia
df_1$exp <- 0.334
df_2$exp <- 1.334
df_3$exp <- 2.334

df_all <- rbind(df_1, df_2, df_3)

df_all <- df_all %>%  
        mutate(distancia = rowid*40)

# eliminar puntos. Luego, puntos con mucha dispercion
# df_all <- df_all %>% 
#     filter(exp == 2.334 & ( distancia == 5480 | distancia == 5520))
# df_all <- df_all %>% 
#     filter(exp == 0.334 & ( distancia == 4640 | distancia == 4680))

# Zr20Nb comienza en 4405 ((127*40)-675)
# Zr comienza en 5755 ((127*40)+675)
g <- df_all %>% 
        pivot_longer(Er:Nb,
                     names_to = "Elementos",
                     values_to = "Intensidad") %>% 
        ggplot(aes(x = distancia, y = Intensidad, colour = Elementos)) +
                geom_point() + geom_line() + 
                # centro del Er 127*40 = 5080um
                geom_vline(xintercept = c((127*40)+675, (127*40)-675)) + 
                annotate(geom = 'text', 
                         x = 5075, y = 0.020, 
                         label = 'Er \n (Inicial)') +
                annotate(geom = 'text', 
                         x = 1800, y = 0.020, 
                         label = 'Zr20Nb \n (Inicial)') +
                annotate(geom = 'text', 
                         x = 7500, y = 0.020, 
                         label = 'Zr \n (Inicial)') +
                xlab('Distancia (micrones)') +
                ylab('Cantidad de cuentas (UA)') +
                theme_test()

g #%>% plotly::ggplotly()

# df_Er_expI <- df_1 %>% select(rowid, Er) %>% mutate(distancia = rowid*40)
# df_Er_expII <- df_2 %>% select(rowid, Er) %>% mutate(distancia = rowid*40)
# df_Er_expIII <- df_3 %>% select(rowid, Er) %>% mutate(distancia = rowid*40)
# 
# write_delim(df_Er_, file = './outputs/dataZr20Nb.txt', delim = " ", col_names = FALSE)

CALCULO DE COEFICIENTE DE ID

df_Er <- df_all %>% 
    select(rowid, Er, distancia) %>% 
    set_names('rowid', 'cuentas', 'distancia')

FUN.perfil <- function(element){
    # element: vector
    element %>% 
        ggplot(aes(x = distancia, y = cuentas)) +
        geom_point() +
        labs(x = 'Distancia (um)', y = 'Numero de cuentas')
}

FUN.perfil(df_Er) %>% plotly::ggplotly()
d_ini <- 6080   # distancia inicial del perfil
df_Er_short <- df_Er %>% 
    filter(distancia >= d_ini & distancia <=7750) %>% 
    mutate(distancia = distancia - d_ini)
    

FUN.perfil(df_Er_short)

erfcinv <- function(y) {
    y[y < 0 | y > 2] <- NA
    -qnorm(y/2)/sqrt(2)
}

#n_max <- max(df_Er_short$cuentas) - 0.0004

# df_Er_short %>% 
#     mutate(erfc_inv = erfcinv(cuentas/n_max)) %>% 
#     ggplot(aes(x = distancia, y = erfc_inv)) +
#     geom_point() +
#     geom_smooth(method = 'lm')

FUN.0intercept <- function(DF, n){
    DF <- DF %>% mutate(erfc_inv = erfcinv(cuentas/n)) 
    fit <- lm(erfc_inv ~ distancia, DF)
    fit$coefficients[[1]]
}

iter <- 1
inter <- 1
paso <- 0.00001
v_int <- numeric()

while(inter > 10^-4){
    if(iter == 1){
        n_max <- max(df_Er_short$cuentas)
        inter <- FUN.0intercept(df_Er_short, n_max)
    }else{
        n_max <- n_max - paso
        inter <- FUN.0intercept(df_Er_short, n_max)
    }
    iter <- iter + 1
    v_int <- c(v_int, inter)
    #if(iter >= 1000) break
}

#plot(v_int)

df_Er_short %>% 
    mutate(erfc_inv = erfcinv(cuentas/n_max)) %>% 
    ggplot(aes(x = distancia, y = erfc_inv)) +
    geom_point() +
    geom_smooth(method = 'lm')
`geom_smooth()` using formula 'y ~ x'

fit <- lm(erfc_inv ~ distancia, df_Er_short %>% 
              mutate(erfc_inv = erfcinv(cuentas/n_max)))
fit

Call:
lm(formula = erfc_inv ~ distancia, data = df_Er_short %>% mutate(erfc_inv = erfcinv(cuentas/n_max)))

Coefficients:
(Intercept)    distancia  
 -0.0014731    0.0001166  
D <- (10^-12)/(4*(1.7021*10^7)*(fit$coefficients[[2]])^2)
D
[1] 1.080918e-12
write_delim(df_Er_short, file = './outputs/dataZr_tresexp.txt', delim = " ", col_names = FALSE)
LS0tDQp0aXRsZTogIlIgTm90ZWJvb2siDQpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sNCi0tLQ0KDQpgYGB7ciBkYXRvc30NCmxpYnJhcnkodGlkeXZlcnNlLCB2ZXJib3NlID0gRkFMU0UpDQoNCkZVTi5yZWFkIDwtIGZ1bmN0aW9uKHBhdGgsIG4pew0KICAgICAgICBNIDwtIG51bWVyaWMoKQ0KICAgICAgICBmb3IoaSBpbiAxOm4pew0KICAgICAgICAgICAgICAgIHNwZWMgPC0gcmVhZF90c3YoZmlsZSA9IHBhc3RlKHBhdGgsICJhIixpLCIub2xzIiwgIHNlcCA9ICcnKSwgDQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBza2lwID0gNiwgc2hvd19jb2xfdHlwZXMgPSBGQUxTRSkgJT4lIA0KICAgICAgICAgICAgICAgICAgICAgICAgc2VsZWN0KENvdW50cykNCiAgICAgICAgICAgICAgICBNIDwtIHJiaW5kKE0sIHNwZWMkQ291bnRzKSANCiAgICAgICAgfQ0KICAgICAgICBNDQp9DQoNCk1fZXNwZWN0cm9zIDwtIHZlY3Rvcihtb2RlID0gJ2xpc3QnLCBsZW5ndGggPSAzKQ0KDQpkaXIgPC0gIi4vZXNwZWN0cm9zL2N1cGxhIFpyLUVyLVpyTmIgcGFyYWxlbGEvRXhwSS8iDQpuX3NwZWMgPC0gbGVuZ3RoKGxpc3QuZmlsZXMoZGlyKSkNCk1fZXNwZWN0cm9zW1sxXV0gPC0gRlVOLnJlYWQoZGlyLCBuX3NwZWMpDQoNCmRpciA8LSAiLi9lc3BlY3Ryb3MvY3VwbGEgWnItRXItWnJOYiBwYXJhbGVsYS9FeHBJSS8iDQpuX3NwZWMgPC0gbGVuZ3RoKGxpc3QuZmlsZXMoZGlyKSkNCk1fZXNwZWN0cm9zW1syXV0gPC0gRlVOLnJlYWQoZGlyLCBuX3NwZWMpDQoNCmRpciA8LSAiLi9lc3BlY3Ryb3MvY3VwbGEgWnItRXItWnJOYiBwYXJhbGVsYS9FeHBJSUkvIg0Kbl9zcGVjIDwtIGxlbmd0aChsaXN0LmZpbGVzKGRpcikpDQpNX2VzcGVjdHJvc1tbM11dIDwtIEZVTi5yZWFkKGRpciwgbl9zcGVjKQ0KDQojIGNhbnRpZGFkIGRlIGFyY2hpdm9zIHBvciBiYXJyaWRvDQpsYXBwbHkoTV9lc3BlY3Ryb3MsIGRpbSkNCg0KIyB3YXZlbGVuZ3RoIHZhbHVlcw0Kd2F2ZWxlbiA8LSByZWFkX3RzdihmaWxlID0gcGFzdGUoZGlyLCAnL2ExLm9scycsIHNlcCA9ICcnKSwgc2tpcCA9IDYsIHNob3dfY29sX3R5cGVzID0gRkFMU0UpICU+JSANCiAgICAgICAgICAgICAgICBzZWxlY3QoV2F2ZWxlbmd0aCkgJT4lIA0KICAgICAgICAgICAgICAgIHJvd2lkX3RvX2NvbHVtbigpDQoNCndhdmVsZW4kV2F2ZWxlbmd0aCA8LSByb3VuZCh3YXZlbGVuJFdhdmVsZW5ndGgsIDQpDQoNCmBgYA0KDQpgYGB7ciBOb3JtYWxpemFjaW9ufQ0KIyBMYSBub3JtYWxpemFjaW9uIGVzIHBvciBkZXRlY3RvciB5IHNlIHJlYWxpemEgcG9yIHN1bWEgdG90YWwuIExvcyBpbmRpY2VzIA0KIyBjb3JyZXNwb25kaWVudGVzIGEgY2FkYSBkZXRlY3RvciBzb246DQojICAgICAgICAgDQojICAgICAgIC0gMToyMDQ4DQojICAgICAgIC0gMjA0OTozOTgzDQojICAgICAgIC0gMzk4NDo1OTI0DQojICAgICAgIC0gNTkyNTo3ODY1DQoNCiMgRnVuY2lvbiBkZSBub3JtYWxpemFjaW9uDQpGVU4ubm9ybSA8LSBmdW5jdGlvbihNLCBuID0gNCwgaW5kZXgpew0KICAgICAgICBsaXN0YSA8LSB2ZWN0b3IobW9kZSA9ICdsaXN0JywgbGVuZ3RoID0gbikNCiAgICAgICAgZm9yKGkgaW4gMTpuKXsNCiAgICAgICAgICAgICAgICAjIHNlcGFyYXIgZGV0ZWN0b3INCiAgICAgICAgICAgICAgICBsaXN0YVtbaV1dIDwtIE1bLGluZGV4W1tpXV1dDQogICAgICAgICAgICAgICAgIyBoYWNlciBtaW5pbW8gPSAwDQogICAgICAgICAgICAgICAgbWluaW1vIDwtIGFwcGx5KGxpc3RhW1tpXV0sIDEsIG1pbikNCiAgICAgICAgICAgICAgICBsaXN0YVtbaV1dIDwtIGxpc3RhW1tpXV0gKyBhYnMobWluaW1vKQ0KICAgICAgICAgICAgICAgICMgRGl2aWRpciBwb3Igc3VtYSB0b3RhbA0KICAgICAgICAgICAgICAgIHRvdGFsIDwtIGFwcGx5KGxpc3RhW1tpXV0sIDEsIHN1bSkNCiAgICAgICAgICAgICAgICBsaXN0YVtbaV1dIDwtIGxpc3RhW1tpXV0gLyB0b3RhbA0KICAgICAgICB9DQogICAgICAgIGxpc3RhICAgICAgIA0KfQ0KDQpudW1fZGV0ZWMgPC0gNCAgICAjIE51bWVybyBkZSBkZXRlY3RvcmVzDQppbmRleF9kZXRlYyA8LSBsaXN0KGMoMToyMDQ4KSwgYygyMDQ5OjM5ODMpLCBjKDM5ODQ6NTkyNCksIGMoNTkyNTo3ODY1KSkgDQoNCk1fbm9ybSA8LSBNX2VzcGVjdHJvcyAlPiUgbWFwKEZVTi5ub3JtLCBuID0gbnVtX2RldGVjLCBpbmRleCA9IGluZGV4X2RldGVjKQ0KDQojIHN1bWEgZWwgZWwgZXNwZWN0cm8gcXVlIGZ1ZSBzZXBhcmFkbyBwYXJhIG5vcm1hbGl6YXIgcG9yIGRldGVjdG9yDQpNX25vcm0gPC0gTV9ub3JtICU+JSBtYXAofiBjYmluZCgueFtbMV1dLCAueFtbMl1dLCAueFtbM11dLCAueFtbNF1dKSkNCg0KbGFwcGx5KE1fbm9ybSwgZGltKSAjIHBhcmEgaW5zcGVjY2lvbg0KDQpgYGANCg0KYGBge3IgZnVuY2lvbmVzfQ0KIyBGdW5jaW9uIHF1ZSBpbnRlZ3JhIGVsIHBpY28NCkZVTi5jdWVudGFzIDwtIGZ1bmN0aW9uKHAsIGJhcnJpZG8sIHNwZWNfaW5pID0gMSwgc3BlY19maW4gPSAyMDApew0KICAgICAgICBpbmQgPC0gd2F2ZWxlbiRyb3dpZFt3aGljaCh3YXZlbGVuJFdhdmVsZW5ndGggPT0gcCldDQogICAgICAgIG0gPC0gTV9ub3JtW1tiYXJyaWRvXV1bc3BlY19pbmk6c3BlY19maW4sKGluZC0yKTooaW5kKzIpXQ0KICAgICAgICBjdWVudGFzIDwtIGFwcGx5KG0sIDEsIHN1bSkNCiAgICAgICAgY3VlbnRhcw0KfQ0KDQojIEZ1bmNpb24gcGFyYSBncmFmaWNhciBwZXJmaWwNCiMgDQpGVU4ucGVyZmlsIDwtIGZ1bmN0aW9uKGRmKXsNCiAgICAgICAgZyA8LSBkZiAlPiUgDQogICAgICAgICAgICAgICAgcGl2b3RfbG9uZ2VyKEVyOk5iLA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICBuYW1lc190byA9ICJFbGVtZW50b3MiLA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICB2YWx1ZXNfdG8gPSAiSW50ZW5zaWRhZCIpICU+JSANCiAgICAgICAgICAgICAgICBnZ3Bsb3QoYWVzKHggPSByb3dpZCo0MCwgeSA9IEludGVuc2lkYWQsIGNvbG91ciA9IEVsZW1lbnRvcykpICsNCiAgICAgICAgICAgICAgICBnZW9tX3BvaW50KCkgIysgZ2VvbV9saW5lKCkNCiAgICAgICAgDQogICAgICAgIGNhcHRpb24gPC0gcGFzdGUoIkxpbmVhcyBlbGVnaWRhczoiLCJcbiIsDQogICAgICAgICAgICAgICAgICAgICAgICAgIkVyOiAiLHBpY29fRXIsIlxuIiwgDQogICAgICAgICAgICAgICAgICAgICAgICAgIlpyOiAiLHBpY29fWnIsIlxuIiwgDQogICAgICAgICAgICAgICAgICAgICAgICAgIk5iOiAiLHBpY29fTmIsIlxuIiwgc2VwID0gJycpDQogICAgICAgIA0KICAgICAgICBnICsgYW5ub3RhdGUoZ2VvbSA9ICd0ZXh0JywgeCA9IDIwLCB5ID0gMC4wMDksIGxhYmVsID0gY2FwdGlvbikgKw0KICAgICAgICAgICAgICAgIHhsYWIoJ0Rpc3RhbmNpYSAobWljcm9uZXMpJykgKw0KICAgICAgICAgICAgICAgIHlsYWIoJ0NhbnRpZGFkIGRlIGN1ZW50YXMgKFVBKScpDQp9DQoNCnBpY29fRXIgPC0gMzY5LjA5MDQNCnBpY29fTmIgPC0gMzIyLjM3OTENCnBpY29fWnIgPC0gMzU3LjA0NDENCmBgYA0KDQpgYGB7ciBQRVJGSUwgQkFSUklETyAxfQ0KIyMjIyMjICMjICMjICBQRVJGSUwgQkFSUklETyAxICAgIyMgIyMgIyMjIyMjDQpiIDwtIDENCmluaSA8LSAxOSAgICNkZXNjYXJ0YSBsb3MgcHJpbWVyb3MgMTggcGVyZmlsZXMNCmZpbiA8LSAyMzQgICN0b21hIGhhc3RhIHBlcmZpbCAyMzQNCkVyIDwtIEZVTi5jdWVudGFzKHBpY29fRXIsIGJhcnJpZG8gPSBiLCBpbmkgLCBmaW4gKQ0KWnIgPC0gRlVOLmN1ZW50YXMocGljb19aciwgYmFycmlkbyA9IGIsIGluaSAsIGZpbiApDQpOYiA8LSBGVU4uY3VlbnRhcyhwaWNvX05iLCBiYXJyaWRvID0gYiwgaW5pICwgZmluICkNCg0KZGZfMSA8LSB0aWJibGUoRXIgPSBFciwgWnIgPSBaciwgTmIgPU5iKSAlPiUgcm93aWRfdG9fY29sdW1uKCkNCmdfYmFycmlkb18xIDwtIEZVTi5wZXJmaWwoZGZfMSkNCmdfYmFycmlkb18xDQpgYGANCg0KYGBge3IgUEVSRklMIEJBUlJJRE8gMn0NCmIgPC0gMg0KaW5pIDwtIDIxDQpmaW4gPC0gbnJvdyhNX25vcm1bW2JdXSkNCkVyIDwtIEZVTi5jdWVudGFzKHBpY29fRXIsIGJhcnJpZG8gPSBiLCBpbmkgLCBmaW4gKQ0KWnIgPC0gRlVOLmN1ZW50YXMocGljb19aciwgYmFycmlkbyA9IGIsIGluaSAsIGZpbiApDQpOYiA8LSBGVU4uY3VlbnRhcyhwaWNvX05iLCBiYXJyaWRvID0gYiwgaW5pICwgZmluICkNCg0KZGZfMiA8LSB0aWJibGUoRXIgPSBFciwgWnIgPSBaciwgTmIgPU5iKSAlPiUgcm93aWRfdG9fY29sdW1uKCkNCmdfYmFycmlkb18yIDwtIEZVTi5wZXJmaWwoZGZfMikNCmdfYmFycmlkb18yDQpgYGANCg0KYGBge3IgUEVSRklMIEJBUlJJRE8gMyB9DQpiIDwtIDMNCmluaSA8LSAyMg0KZmluIDwtIG5yb3coTV9ub3JtW1tiXV0pDQpFciA8LSBGVU4uY3VlbnRhcyhwaWNvX0VyLCBiYXJyaWRvID0gYiwgaW5pICwgZmluICkNClpyIDwtIEZVTi5jdWVudGFzKHBpY29fWnIsIGJhcnJpZG8gPSBiLCBpbmkgLCBmaW4gKQ0KTmIgPC0gRlVOLmN1ZW50YXMocGljb19OYiwgYmFycmlkbyA9IGIsIGluaSAsIGZpbiApDQoNCmRmXzMgPC0gdGliYmxlKEVyID0gRXIsIFpyID0gWnIsIE5iID1OYikgJT4lIHJvd2lkX3RvX2NvbHVtbigpDQpnX2JhcnJpZG9fMyA8LSBGVU4ucGVyZmlsKGRmXzMpDQpnX2JhcnJpZG9fMw0KYGBgDQoNCsK/Q3VhbCBzZXLDoSBsYSByZWZlcmVuY2lhIHBhcmEgdW5pciBsb3MgZGF0b3M/DQoNCi0gICBFbCBFciB0aWVuZSBtdWNoYSBkaXNwZXJjaW9uDQoNCi0gICBUb23DqSBkZSByZWZlcmVuY2lhIGFsIFpyDQoNCmBgYHtyIFN1bWFyIGxvcyB0cmVzIGJhcnJpZG9zfQ0KIyBkZXNwbGF6YW5kbyBpbmRpY2VzDQojIHJvd2lkIHZpZW5lIGEgc2VyIGVsIGdpcm8gZGVsIG1pY3JvbWV0cm8uIDQwIG1pY3JvbmVzDQpkZl8yJHJvd2lkIDwtIGRmXzIkcm93aWQgKyAzDQpkZl8zJHJvd2lkIDwtIGRmXzMkcm93aWQgKyA2DQoNCiMgZXN0byBoYXkgcXVlIGNhbWJpYXJsbyBhIGRpc3RhbmNpYQ0KZGZfMSRleHAgPC0gMC4zMzQNCmRmXzIkZXhwIDwtIDEuMzM0DQpkZl8zJGV4cCA8LSAyLjMzNA0KDQpkZl9hbGwgPC0gcmJpbmQoZGZfMSwgZGZfMiwgZGZfMykNCg0KZGZfYWxsIDwtIGRmX2FsbCAlPiUgIA0KICAgICAgICBtdXRhdGUoZGlzdGFuY2lhID0gcm93aWQqNDApDQoNCiMgZWxpbWluYXIgcHVudG9zLiBMdWVnbywgcHVudG9zIGNvbiBtdWNoYSBkaXNwZXJjaW9uDQojIGRmX2FsbCA8LSBkZl9hbGwgJT4lIA0KIyAgICAgZmlsdGVyKGV4cCA9PSAyLjMzNCAmICggZGlzdGFuY2lhID09IDU0ODAgfCBkaXN0YW5jaWEgPT0gNTUyMCkpDQojIGRmX2FsbCA8LSBkZl9hbGwgJT4lIA0KIyAgICAgZmlsdGVyKGV4cCA9PSAwLjMzNCAmICggZGlzdGFuY2lhID09IDQ2NDAgfCBkaXN0YW5jaWEgPT0gNDY4MCkpDQoNCiMgWnIyME5iIGNvbWllbnphIGVuIDQ0MDUgKCgxMjcqNDApLTY3NSkNCiMgWnIgY29taWVuemEgZW4gNTc1NSAoKDEyNyo0MCkrNjc1KQ0KZyA8LSBkZl9hbGwgJT4lIA0KICAgICAgICBwaXZvdF9sb25nZXIoRXI6TmIsDQogICAgICAgICAgICAgICAgICAgICBuYW1lc190byA9ICJFbGVtZW50b3MiLA0KICAgICAgICAgICAgICAgICAgICAgdmFsdWVzX3RvID0gIkludGVuc2lkYWQiKSAlPiUgDQogICAgICAgIGdncGxvdChhZXMoeCA9IGRpc3RhbmNpYSwgeSA9IEludGVuc2lkYWQsIGNvbG91ciA9IEVsZW1lbnRvcykpICsNCiAgICAgICAgICAgICAgICBnZW9tX3BvaW50KCkgKyBnZW9tX2xpbmUoKSArIA0KICAgICAgICAgICAgICAgICMgY2VudHJvIGRlbCBFciAxMjcqNDAgPSA1MDgwdW0NCiAgICAgICAgICAgICAgICBnZW9tX3ZsaW5lKHhpbnRlcmNlcHQgPSBjKCgxMjcqNDApKzY3NSwgKDEyNyo0MCktNjc1KSkgKyANCiAgICAgICAgICAgICAgICBhbm5vdGF0ZShnZW9tID0gJ3RleHQnLCANCiAgICAgICAgICAgICAgICAgICAgICAgICB4ID0gNTA3NSwgeSA9IDAuMDIwLCANCiAgICAgICAgICAgICAgICAgICAgICAgICBsYWJlbCA9ICdFciBcbiAoSW5pY2lhbCknKSArDQogICAgICAgICAgICAgICAgYW5ub3RhdGUoZ2VvbSA9ICd0ZXh0JywgDQogICAgICAgICAgICAgICAgICAgICAgICAgeCA9IDE4MDAsIHkgPSAwLjAyMCwgDQogICAgICAgICAgICAgICAgICAgICAgICAgbGFiZWwgPSAnWnIyME5iIFxuIChJbmljaWFsKScpICsNCiAgICAgICAgICAgICAgICBhbm5vdGF0ZShnZW9tID0gJ3RleHQnLCANCiAgICAgICAgICAgICAgICAgICAgICAgICB4ID0gNzUwMCwgeSA9IDAuMDIwLCANCiAgICAgICAgICAgICAgICAgICAgICAgICBsYWJlbCA9ICdaciBcbiAoSW5pY2lhbCknKSArDQogICAgICAgICAgICAgICAgeGxhYignRGlzdGFuY2lhIChtaWNyb25lcyknKSArDQogICAgICAgICAgICAgICAgeWxhYignQ2FudGlkYWQgZGUgY3VlbnRhcyAoVUEpJykgKw0KICAgICAgICAgICAgICAgIHRoZW1lX3Rlc3QoKQ0KDQpnICMlPiUgcGxvdGx5OjpnZ3Bsb3RseSgpDQpgYGANCg0KYGBge3Igc2F2ZSBkYXRhfQ0KIyBkZl9Fcl9leHBJIDwtIGRmXzEgJT4lIHNlbGVjdChyb3dpZCwgRXIpICU+JSBtdXRhdGUoZGlzdGFuY2lhID0gcm93aWQqNDApDQojIGRmX0VyX2V4cElJIDwtIGRmXzIgJT4lIHNlbGVjdChyb3dpZCwgRXIpICU+JSBtdXRhdGUoZGlzdGFuY2lhID0gcm93aWQqNDApDQojIGRmX0VyX2V4cElJSSA8LSBkZl8zICU+JSBzZWxlY3Qocm93aWQsIEVyKSAlPiUgbXV0YXRlKGRpc3RhbmNpYSA9IHJvd2lkKjQwKQ0KIyANCiMgd3JpdGVfZGVsaW0oZGZfRXJfLCBmaWxlID0gJy4vb3V0cHV0cy9kYXRhWnIyME5iLnR4dCcsIGRlbGltID0gIiAiLCBjb2xfbmFtZXMgPSBGQUxTRSkNCmBgYA0KDQojIyMgQ0FMQ1VMTyBERSBDT0VGSUNJRU5URSBERSBJRA0KDQpgYGB7cn0NCmRmX0VyIDwtIGRmX2FsbCAlPiUgDQogICAgc2VsZWN0KHJvd2lkLCBFciwgZGlzdGFuY2lhKSAlPiUgDQogICAgc2V0X25hbWVzKCdyb3dpZCcsICdjdWVudGFzJywgJ2Rpc3RhbmNpYScpDQoNCkZVTi5wZXJmaWwgPC0gZnVuY3Rpb24oZWxlbWVudCl7DQogICAgIyBlbGVtZW50OiB2ZWN0b3INCiAgICBlbGVtZW50ICU+JSANCiAgICAgICAgZ2dwbG90KGFlcyh4ID0gZGlzdGFuY2lhLCB5ID0gY3VlbnRhcykpICsNCiAgICAgICAgZ2VvbV9wb2ludCgpICsNCiAgICAgICAgbGFicyh4ID0gJ0Rpc3RhbmNpYSAodW0pJywgeSA9ICdOdW1lcm8gZGUgY3VlbnRhcycpDQp9DQoNCkZVTi5wZXJmaWwoZGZfRXIpICU+JSBwbG90bHk6OmdncGxvdGx5KCkNCmBgYA0KDQpgYGB7cn0NCmRfaW5pIDwtIDYwODAgICAjIGRpc3RhbmNpYSBpbmljaWFsIGRlbCBwZXJmaWwNCmRmX0VyX3Nob3J0IDwtIGRmX0VyICU+JSANCiAgICBmaWx0ZXIoZGlzdGFuY2lhID49IGRfaW5pICYgZGlzdGFuY2lhIDw9Nzc1MCkgJT4lIA0KICAgIG11dGF0ZShkaXN0YW5jaWEgPSBkaXN0YW5jaWEgLSBkX2luaSkNCiAgICANCg0KRlVOLnBlcmZpbChkZl9Fcl9zaG9ydCkNCmBgYA0KDQpgYGB7cn0NCmVyZmNpbnYgPC0gZnVuY3Rpb24oeSkgew0KICAgIHlbeSA8IDAgfCB5ID4gMl0gPC0gTkENCiAgICAtcW5vcm0oeS8yKS9zcXJ0KDIpDQp9DQoNCiNuX21heCA8LSBtYXgoZGZfRXJfc2hvcnQkY3VlbnRhcykgLSAwLjAwMDQNCg0KIyBkZl9Fcl9zaG9ydCAlPiUgDQojICAgICBtdXRhdGUoZXJmY19pbnYgPSBlcmZjaW52KGN1ZW50YXMvbl9tYXgpKSAlPiUgDQojICAgICBnZ3Bsb3QoYWVzKHggPSBkaXN0YW5jaWEsIHkgPSBlcmZjX2ludikpICsNCiMgICAgIGdlb21fcG9pbnQoKSArDQojICAgICBnZW9tX3Ntb290aChtZXRob2QgPSAnbG0nKQ0KDQpGVU4uMGludGVyY2VwdCA8LSBmdW5jdGlvbihERiwgbil7DQogICAgREYgPC0gREYgJT4lIG11dGF0ZShlcmZjX2ludiA9IGVyZmNpbnYoY3VlbnRhcy9uKSkgDQogICAgZml0IDwtIGxtKGVyZmNfaW52IH4gZGlzdGFuY2lhLCBERikNCiAgICBmaXQkY29lZmZpY2llbnRzW1sxXV0NCn0NCg0KaXRlciA8LSAxDQppbnRlciA8LSAxDQpwYXNvIDwtIDAuMDAwMDENCnZfaW50IDwtIG51bWVyaWMoKQ0KDQp3aGlsZShpbnRlciA+IDEwXi00KXsNCiAgICBpZihpdGVyID09IDEpew0KICAgICAgICBuX21heCA8LSBtYXgoZGZfRXJfc2hvcnQkY3VlbnRhcykNCiAgICAgICAgaW50ZXIgPC0gRlVOLjBpbnRlcmNlcHQoZGZfRXJfc2hvcnQsIG5fbWF4KQ0KICAgIH1lbHNlew0KICAgICAgICBuX21heCA8LSBuX21heCAtIHBhc28NCiAgICAgICAgaW50ZXIgPC0gRlVOLjBpbnRlcmNlcHQoZGZfRXJfc2hvcnQsIG5fbWF4KQ0KICAgIH0NCiAgICBpdGVyIDwtIGl0ZXIgKyAxDQogICAgdl9pbnQgPC0gYyh2X2ludCwgaW50ZXIpDQogICAgI2lmKGl0ZXIgPj0gMTAwMCkgYnJlYWsNCn0NCg0KI3Bsb3Qodl9pbnQpDQoNCmRmX0VyX3Nob3J0ICU+JSANCiAgICBtdXRhdGUoZXJmY19pbnYgPSBlcmZjaW52KGN1ZW50YXMvbl9tYXgpKSAlPiUgDQogICAgZ2dwbG90KGFlcyh4ID0gZGlzdGFuY2lhLCB5ID0gZXJmY19pbnYpKSArDQogICAgZ2VvbV9wb2ludCgpICsNCiAgICBnZW9tX3Ntb290aChtZXRob2QgPSAnbG0nKQ0KYGBgDQoNCmBgYHtyfQ0KZml0IDwtIGxtKGVyZmNfaW52IH4gZGlzdGFuY2lhLCBkZl9Fcl9zaG9ydCAlPiUgDQogICAgICAgICAgICAgIG11dGF0ZShlcmZjX2ludiA9IGVyZmNpbnYoY3VlbnRhcy9uX21heCkpKQ0KZml0DQpgYGANCg0KYGBge3J9DQpEIDwtICgxMF4tMTIpLyg0KigxLjcwMjEqMTBeNykqKGZpdCRjb2VmZmljaWVudHNbWzJdXSleMikNCkQNCmBgYA0KDQpgYGB7cn0NCndyaXRlX2RlbGltKGRmX0VyX3Nob3J0LCBmaWxlID0gJy4vb3V0cHV0cy9kYXRhWnJfdHJlc2V4cC50eHQnLCBkZWxpbSA9ICIgIiwgY29sX25hbWVzID0gRkFMU0UpDQpgYGANCg==